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Abstract 


A numerical method for calculating unsteady two-dimensional boundary layers in 
incompressible laminar and turbulent flows is described, and applied to a single airfoil 
changing its incidence angle in time. The solution procedure adopts a first order panel 
method with a simple wake model to solve for the inviscid part of the flow', and an implicit 
finite difference method for the viscous part of the flow. Both procedures integrate in time 
in a step-by-step fashion, in course of which each step involves the solution of the elliptic 
Laplace equation and the solution of the parabolic boundary layer equations. The Reynolds 
shear stress term of the boundary layer equations is modeled by an algebraic eddy viscosity 
closure, which essentially is the steady Cebeci-Smith model. The location of transition is 
predicted by an empirical data correlation originating from Michel. Since transition and 
turbulence modeling are key factors in the prediction of viscous Tow's, their accuracy will 
be of dominant influence to the overall results. 

The presented methodology offers a decisive cost advantage as compared to solvers of the 
full Navier-Stokes equations. The efficiency of our method derives from its viscous and 
inviscid components, which both bypass the far more expensive field solution. Conceptua y 
the method makes no assumptions as to the type of motion or the shape of the airfoil, and 
is subject only to the limitations imposed by the boundary layer approximation and by the 
assumption of an irrotational inviscid flow field. As a consequence of the currently em- 
ploved boundarv layer approach, in which the pressure is prescribed, the method cannot 
cope with the problem of flow separation. Future work is planned to address the problem 
of unsteady flow separation, W'hich will require an interactive coupling of the inviscid and 
viscous flow solvers. 
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1. Introduction 

The prediction of unsteady viscous flow phenomena is of importance for a variety of aero- 
dynamic devices, including turbomachinery bladings, helicopter rotor blades, wind energy 
devices, and rapidly maneuvering aircraft. The unsteadiness of the problem can result either 
from an unsteady onset flow acting on a stationary body, or from an unsteadily moving 
body in an initially uniform onset flow. Examples for the first type of problems are the flow 
through the stator of a turbomachine being exposed to the w r akes of the previous stages and 
the flow past aircraft travelling through atmospheric gusts. Examples for the latter type are 
the flows past rapidly maneuvering aircraft in quiet air and the flow through the vibrating 
blades of a turbomachine. Manv problems of practical interest involve both forms of un- 
steadiness, like the flowfield surrounding a helicopter rotor, whose blades change both 
speed and incidence angle and are exposed to the wakes of the previous blades. Most of 
todays numerical tools hardly can keep up with the complexity of these flows. The majority 
of today's design methods either assumes inviscid flow conditions or relies upon empirical 
data correlations. However, the advances in computer technology have initiated the devel- 
opment of methods, which do take account of viscous effects. Large efforts are now under 
way in order to attempt the solution of the full Navier-Stokes equations or of some 
asymptotic approximation to the Navier-Stokes equations. 

Although the full Navier-Stokes equations are generally accepted as the exact physical 
model of the flow, asymptotic approximations continue to be pursued because of their rel- 
ative merits with regard to computer time and storage and because of their ability of pro- 
moting the understanding of flows. Of course, asymptotic approximations cannot be 
expected to be universally applicable and their inadequate use can lead to numerical 
breakdowns. The Goldstein singularity /22/ in steady boundary' layers and Van Dommelen s 
singularity /43 / in unsteady boundary layers are prominent examples of this very fact. Yet 
the appearance of singularities can be considered to be an advantage, since they often reveal 
the transition to another flow regime and indicate the need for another asymptotic de- 
scription. Under a disadvantageous scenario the solution of an asymptotic approximation 
mieht- simply diverge from the Navier-Stokes solution without giving any indication. In 
spite of the inherent limitations of approximations it is the simplified model that allows us 
to sain fundamental knowledge by reducing the problem to its essential. 

The Navier-Stokes equations are of the elliptic type, which means that information propa- 
gates in all directions. In principle, information is transmitted via pressure, via the viscous 
(and turbulent) stress gradients, and via the local velocity. While the signals of the first tw o 
types spread in all directions, the signals of the third are convected with the stream. A sol- 
ution of the full Navier-Stokes equations requires thus a simultaneous processing of the 
pressure, the velocity components, and the stress tensor throughout the flowfield, which for 
the present rules out this procedure as an engineering tool because of costs. Simplifications 
of the Navier-Stokes equations consist of restricting the signal propagation in certain parts 
of the flowfield. Provided that the longitudinal stress gradients are small compared to the 
cross-stream shear gradients, the signals due to the viscous stresses can be assumed to not 
propagate in the streamw'ise direction. In practice, this is achieved by omitting these minor 
stress terms, which leads to the so-called parabolized or thin-layer Navier-Stokes equations. 
The first name refers to the partially changed characteristics of the new’ equations, although 
the elliptic behavior is retained due" to the pressure signals being sent in all directions. The 
second name indicates that the viscous effects ought to be confined to a wall-bounded re- 
gion, whose cross-stream dimension is small. In contrast to the boundary layer equations 
this approach retains the momentum equation across the viscous layer, so that the strong 
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coupling between inner viscous and outer inviscid flow is maintained. The classical bound- 
ary layer approximation goes one step further: In addition to limiting the viscous effects to 
a thin layer, the calculation of the pressure is separated from the calculation of the viscous 
flowfield. This leaves one flowfield, in which there are no signals due to viscous stress gra- 
dients, and another one, in which there are no pressure signals. The flowfield that deter- 
mines the pressure can be described by the Euler equations, which evolve from the 
Navier-Stokes equations by omitting all viscous terms. The viscous layer is governed by 
Prandtl's boundary layer equations, in which longitudinal diffusion' terms have been 
dropped and which require the pressure being imposed from outside. This simplest set of 
equations can be developed under the assumption that inside the boundary layer the gra- 
dient of any velocity component in the cross-stream direction is much larger than its gra- 
dient in the streamwise direction. It is this very condition, which makes the classical 
boundary layer equations so effective and at the same time so vulnerable. Because the sol- 
ution of the boundary layer equations can be obtained by progressing step by step in the 
streamwise direction (a consequence of the parabolic characteristics), boundary layers are 
by an order of magnitude less costly to solve than viscous llowflelds governed' by the 
Navier-Stokes equations. Because the outer flowfield has been decoupled from the viscous 
layer (a consequence of the dropped momentum equation in the cross-stream direction and 
of the imposed pressure on the boundary layer), the boundary layer approach does not 
apply in situations where the viscous layer has a substantial effect on the outer inviscid 
flowfield. This prevents regions of separated flow from successful treatment by means of the 
classical boundary layer theory. 

Since the appearance of the boundary layer theory in 1904, many attempts have been 
undertaken to turn the boundary layer approach into an even more versatile tool and to 
remove some of its limitations. It was found that the Goldstein singularity, which occurs 
at the point of zero wall shear in steady boundary layers subject to a prescribed pressure 
distribution, can be removed by prescribing instead the displacement thickness (or the wall 
shear / 5/). While in direct methods the pressure is determined by the inviscid flow and the 
displacement thickness by the viscous flow, such a boundary condition would swap the 
roles of pressure and displacement thickness. Thus inverse methods determine the velocity 
at the boundary layer edge, and consequently the pressure along the boundary layer, as part 
of the viscous flow solution. But inverse methods introduced new problems (slow conver- 
gence), and also common sense suggests that the flow over an airfoil involves both types 
ol interaction. This led to the development of the interactive boundary layer methods, 
which are designed to ensure interaction in both ways, that is the viscous flow is allowed 
to influence the inviscid flow and vice versa. The essence of these methods is an interactive 
boundary condition, which prescribes at the boundary layer edge a linear combination of 
the unknown velocity and displacement thickness. The overall solution is now obtained in 
an iterative process, which accounts for the displacement effect upon the outer inviscid 
flow. So far, interactive boundary layer procedures have an outstanding record in dealing 
with such diversified problems like subsonic and transonic flows /29/, single airfoil /IO /, 
multi-element airfoil and airfoil-spoiler flows / 29/, steady /45/ and unsteady flows 120,22 f. 
Interactive boundary layer procedures can be considered as inexpensive alternatives to the 
Navier-Stokes solvers, since they are, as well, capable of handling (moderatelv) separated 
flows. 

Concepts similar to those of the steady boundary layer can be applied to the unsteady 
problem. Unfortunately, our knowledge of unsteady boundary layers is less sound and the 
status of unsteady procedures is less settled than that of steady ones. For example, the cri- 
terion of unsteady flow separation is still a matter of discussion. While in steady boundary 
layers the point of separation coincides with the point of zero wall shear, such a simple 
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correlation does not exist for unsteady boundary layers. There have been other separation 
criteria proposed (for a review see /4l / or /47 /), like the Moore-Rott-Sears criterion, which 
derives from a boundary layer along a moving wall, but their universal applicability remains 
disputed. It is, however, clear that the point of zero wall shear, in general, will be nonsin- 
gular in unsteady boundary layers, and consequently direct methods do not experience a 
Goldstein-like singularity at this very point. But unsteady boundary layers are by no means 
free of singularities: With the external velocity as prescribed datum the unsteady boundary 
layer of an impulsively started cylinder develops the so-called Van Dommelen singularity 
(which is not of the Goldstein type) within a finite time. Since the discovery of this 
singularity (in 1977) numerous attempts have been made to remove this singularity, but 
none of them could achieve more than to delay its appearance. Certainly, the development 
of unsteady interactive boundary layer procedures is only at its beginning, and further work 
is needed in order to better understand unsteady flow separation. 

As yet another problem of great importance, concerning both Navier-Stokes methods and 
boundary layer methods alike, has not been solved in a satisfactory way. Both Navier- 
Stokes equations and asymptotic approximations do not contain enough information about 
turbulence to form a closed set of equations. The unknown Reynolds stresses must be 
modeled empirically, so that the number of unknowns is reduced to equal the number of 
equations. The most common approach is to define an eddy viscosity, which may be ob- 
tained from the solution of either algebraic equations or differential equations, like the 
transport equation of turbulence energy. Algebraic models offer simplicity and accuracy 
over the range of data upon which they are based, while transport models are usually ap- 
plicable to a wider range of flows, such as flows with a strong adverse pressure gradient and 
extensively separated flow's. But not only the modeling of turbulence, also the prediction 
of transition holds plenty of uncertainties. Since flows past airfoils, in particular at low 
Reynolds numbers, comprise laminar, transitional, and turbulent regimes, there is the need 
to specify both the onset of transition and the extent of the region between laminar and 
turbulent flow. The choice lies here between empirical data correlations, which relate the 
locus of transition to some mean flow' parameter, and linear stability theory, which requires 
the solution of the linear stability equations and subsequently correlates the onset of tran- 
sition to the amplification ratio of disturbances. From a practical point of view the former 
methods are advantageous because of the ease of numerical implementation, w'hile the lat- 
ter ones offer a greater degree of generality. Since both transition and turbulence have a 
dominant effect on the flow results, their importance hardly can be overestimated and their 
underlying assumptions cannot be scrutinized enough. A further approximation of all en- 
gineering methods requires that the time-accurate equations are always replaced by time- 
averaged equations. For the time being the solution of the full Navier-Stokes equations 
with a time-scale that would resolve the motion of turbulent eddies is simply not feasible 
with a reasonable employment of computer resources. Therefore the Navier-Stokes 
equations and all the asymptotic approximations are time-averaged over a time-scale that 
is too larse to resolve any scales of the motion in turbulent flow', yet is small compared to 
the aerodynamic time-scale of the unsteady flow itself. So if we, earlier in this chapter, re- 
ferred to procedures to solve the Navier-Stokes equations or some reduced form of them, 
we implicitly meant the solution of the Reynolds-averaged equations, which do not resolve 
turbulence. 

Finally, w'e want to focus the discussion on the approach pursued here. The subject, we 
examine, is the unsteady two-dimensional flow of an incompressible viscous fluid past a 
single airfoil. The fluid is assumed to be initially at rest or in uniform motion and the un- 
steadiness of the flow shall result solely from the unsteadily moving airfoil. There are no 
assumptions made as to the type of motion the airfoil is executing or as to the shape of the 
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airfoil. The chosen methodology attempts the solution of the complete flow by a synthesis 
of an outer potential flow and an inner boundary layer flow. The solution "of the outer 
inviscid flowfield is obtained by an unsteady panel method, which, in contrast to steady 
methods, has to model the continuous shedding of vorticity into the trailing wake. As in the 
original panel method by Smith and Hess /24/, the singularity distribution is made up by 
uniform source panels and a uniform vorticity around the airfoil. In addition, the wake 
consists of a vorticity panel, which is attached to the trailing edge, and a series of point 
vortices, which are convected downstream with the local particle speed. The solution of the 
viscous flow is determined by a finite difference approximation to the boundary layer 
equations cast in dimensionless, physically transformed variables. In order to properly ac- 
count for flow reversals, which in unsteady flows are not necessarily related to flow sepa- 
rations, the currently employed method utilizes a numerical technique advocated by Cebeci, 
that is the use of both the regular box scheme /9/ and the characteristic box scheme 
/7, 12, 1 3/. In addition, the solution of the unsteady boundary layer equations requires that 
initial and upstream boundary conditions be provided. On the assumption that the time 
derivatives of all variables vanish at the initial time (of the computational domain), the in- 
itial conditions can be obtained from a steady-state solution. The upstream boundary con- 
ditions are generated by a procedure, involving the iterative use of the characteristic box 
scheme, which allows for both the movement of the stagnation point and the eventual oc- 
currence of flow reversals. The Reynolds shear stress term is approximated by means of the 
eddy viscosity concept, using the algebraic formulation of the Cebeci-Smith model / 6/. The 
onset of transition is predicted by an empirical data correlation originating from Michel 
1331 , while the finite region of transition is taken into account by an intermittency distrib- 
ution stemming from Chen and Thyson / 1 4 /. As of now the method performs calculations 
in the direct mode only, that is the outer boundary condition prescribes the external veloc- 
ity, which, prevents the method from successfully dealing with separated flows. The inter- 
active coupling of the potential flow and the boundary layer flow, which shall permit the 
inclusion of separated flow regions, is going to be addressed in the near future. 

The present report is divided into five chapters, the first of which is this introduction. The 
second chapter is concerned with the formulation of the unsteady panel method, in partic- 
ular, attention is given to the unsteady pressure equation, the unsteady Kutta condition, 
and the wake model. The third chapter is devoted to the description of the unsteady 
boundary layer method. In order to gain a better understanding of the numerical necessi- 
ties, first the mechanism of signal propagation in unsteady boundary layers is outlined. 
Subsequently, the finite difference approximation to the unsteady boundary layer equations 
is introduced, and the algorithms of the numerical solution are discussed "briefly. A rather 
short account is then given of the modeling of turbulence and of the prediction of transi- 
tion. A preliminary reflection on unsteady flow separation concludes the third chapter. The 
results of the flow past a single airfoil undergoing a ramp-change in angle of attack are 
presented in the fourth chapter. The final and fifth chapter identifies summary conclusions 
and addresses some of the goals of future work. The reader is alerted to the fact that the 
present report is of an interim nature covering the period from September 1986 through 
August 1987. 
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2. The unsteady potential flow method 

The calculation of the inviscid incompressible flow over an airfoil, which is unde^oing an 
arbitrary time-dependent motion, is performed by a panel method that includes a model for 
the continuous shedding of vorticity into the trading wake. The existence of a ' sh * 
behind the airfoil can be explained by the Helmholtz theorem, which requires that an 
change in the circulation around the airfod must be matched by the appearance of an equal 
counter^ tex om Xre in the flowfield. Naturally, the appearance of this countervortex 
takes place at the trading edge of the airfoil. Since free vortices move with the particles of 
the surrounding fluid, all countervortices together will form a sheet, which extends Tronr the 
trailing edge in the downstream direction. Also the bound vortex ol an airfoil in - 
is balanced by a series of countervortices, which has been generated when the airfoil started 
its motion. But there the existence of these countervortices does not present any difficulties 
because the countervortices had been carried too far downstream to have any effect on the 
flow at^ the aS Nofso in an unsteady flow, in which the airfod may experience a con- 
tinuous change of circulation: The countervortices will stay for some time close enou e h to 
affect the flow over the airfoil, though their influence will lessen as they travel downstream. 
The flow over the airfoil in turn will influence the process of vortex shedding and the shape 
of theTortlx sheet in the near wake, while the flow in the far wake is determined maiffiy 
bv the vorticitv which is stored up there. The presence of the countervortices provides th 
flow with kind of a memory in that the flow at a particular time is affected by the boun 
circulation of the past. Anv mathematical model that describes unsteady flows must hence 
be equipped with a mechanism to simulate the process of shedding vorticity into the wak , 
which distinguishes the unsteady flow from its steady counterpart. 

The inviscid flow solution is based on four principles, namely the conservation of mass, the 
Conservation of momentum, the conservation of 

flow. The basic governing equation, comprising the contmuity equation 
incompressible fluid (div V = 0) and the condition of irrotational flow (curl V = 0), is the 
Laplace equation 

div (grad (p) = 0 

with (p(X,Y, t) denoting the potential for the velocity. Irrotationality is the necessary and 

sufficient condition for the existence of a velocity potential such that V — grad <p. 
Interestingly the Laplace equation applies to steady and unsteady flows alike, 
which unsteady flow methods can be obtained by extending steady 
significant feature of the Laplace equation is that the velocity turns out to be decoupiea 
from 1 the ^pressure and hence the calculation of velocity and pressure can be performed 
separatelv P and consecutively. Following the calculation of the velocity flowfield, the ^pres- 
sure can be determined from the Bernoulli equation, which derives from the law of conser- 
vation of momentum under the previously stated assumptions 

( 2 . 2 ) 


P~ Pc 


iL+M- — o 

2 dt 


This form of the Bernoulli equation refers to an inertial frame of reference \_X,Y], in . w 1 
the motion of the fluid is described by the velocity V(X,Y, /). Since it is more convemen 
pose the problem in a frame of reference that is attached to the moving airfoil, the Bernoulli 

eauation must be written for a moving frame of reference. With K s (x,y, r) denoting the 
velocity of a point in the moving system [*,y] with respect to the inertial system, the ve- 
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locity ujx,y, t) of a fluid particle as viewed in the moving frame can be related to the ve- 
locity V{x,y, i) of the same fluid panicle as viewed in the inertial frame according to 


y — Vs + 0 ~ ^ & (xj —y i ) + o 


( 2 . 3 ) 


The velocity fs(jr,y, 0 can funher be decomposed into a translational speed V ( r ) and a 

TJT a l SV ? d ’ Q{[) • the f0ll0win g 11 4 be a ^™d that ^t^traSlaflo^ 
includes the free stream , so that the fluid panicles are initially at rest with respect to 

nh e JHt rtl fh SyStem ‘u S thC 3irf °i a PP roaches a particular station in the inertial system an 
observer there could perceive the initially motionless fluid particles making wav for the 

mIfim? OUgh th e fluid particles might neither return to their original position nor be- 
come stationary- within a finite time, they do not leave their original neighborhood. Because 
fluid particles (as viewed in the inertial system) do only get disturbed from their original 
position, V and 4 » will be called the disturbance velocity and the disturbance potential re- 
spectively. But while the velocity can be given with respect to any frame of reference this 
is not the case for the potential. Vorticity and circulation are, as well as the velocity ’con- 
cepts relative to a frame of reference, and therefore the same flow, viewed from different 
frames, may well be either irrotational or rotational. This, in general, will not allow to in- 
troduce a velocity potential for the flowfield as seen in the moving system. But we have the 
option to write the velocity potential in terms of other independent variables. So if we re- 
place the disturbance potential <f>(X, Y,t) by the disturbance potential <p(x,y, t) , which de- 

moving frame ^ n ° Wfield ^^Sradq>), we arrive at the Bernoulli equation for the 


P-Poo Vs 0 dq> 

P 2 + 2 + dt 


= 0 


( 2 . 4 ) 

fin H !n wm , M e S Uat1 ^ and the , Bernoum ec f uatl0n ( with a zero time-derivative of the po- 
tential) w ould be sufficient to determine steady flow. Unsteady flow calls for one additional 

pinciple, the conservation of vorticity. The Helmholtz theorem requires that the totafdS. 

^ 10 V n a potential flowfield must be preserved, which, interpreted for an airfoil flow 

the ch r nge m . circulatI0n around the airfoil must be matched bv an equal and 
opposite change of vorticity m the wake. ' 4 “^ <mu 

IS^n C r nSider the , b< ?“ ndar y conditions, which complement the governing equations in 
order to form a solvable problem. There are two of them, namelv the condiS nf 
tangential flow and the Kutta condition. The condition of tangential flow is easilv stated 
° f Re moving frame of reference, that is the velocity o of a fluid particle with re- 
Sce°of ‘,h“woU mOVm8 ° f reference) have a componenfnorL ,o *e 


i> • n = 0 


( 2 . 5 ) 


n addition, another boundary condition is needed to uniquely establish a solution It is 
common practice to make use of the empirical Kutta condition for that purpose. The Kutta 
con ltion postulates zero load at the trailing edge in order to ensure that the flows of upper 

Tnv iTn r SUrfaC K the airfoi i srnoothi .y- The trouble with this procedure is that an 
• 1 met hod makes use of a condition that derives from an experimental observa- 

tion in a viscous fluid. I t is thus the consideration of viscosity, which allows to fix the cir- 
culation around the airfoil and consequently yields a unique solution of the inviscid flow 
pro lem. While in steady flows the use of the Kutta condition is well accepted there is ex- 
perimental evidence that the condition of zero load at the trailing edge is violated in certain 
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unsteady flows. But in the absence of anything better most - if not all - unsteady methods 
(includin'* the present approach) stick to the zero load trailing edge condition. In gen 
such a condition will not lead to a stagnation point at the trailing edge nor 'to [ velocitl ^ s 
of equal magnitude on upper and lower surface at the trailing edge. Only when the rate ot 
change in circulation is zero, that is steady flow conditions, will the velocities approach the 
same" value. But this is, on physical grounds, more acceptable than having a nonzero pres- 
sure differential at the trailing edge, and furthermore the zero load condition does not 
conflict with the shedding of vorticity at the trailing edge. 

We will now turn to the numerical technique to solve the unsteady potential now B ec ause 
the principle of superposition applies to the linear Laplace equation, the flowfield can be 
builf up bv simple elementary flows. The significance of elementary solutions lies in that 
they automatically satisfy the Laplace equation and the task left is to compose them such 
that thev also satisfy the boundary conditions. For obvious reasons numerical techniq 
seek to utilize the simplest elementary flows available, which are those of a point source, a 
point sink, and a point vortex. There are many legitimate choices of these elementary sol- 
utions (also termed singularities), and there are a great many ways to distribute them, 
pmferabfv however, on or close to the surface of the airfoil. Correspondingly, numerous 
Techniques are now available to achieve the numerical solution of an m^iscid incompressible 
flow. The most popular technique is associated with Smith and Hess /24/, who devised he 
so-called panel method in the early sixties. The first extension of the panel method to 
unsteady problem was given by Giesing in 1968 /21/, and since then a coup e of authors, 
amongst them Basu and Hancock /2/, and more recently Kim and Mook / 27 / followed 
proposing some variations of the unsteady approach. With the exception of the wake, 
which needs to be simulated in the unsteady flow model, steady and unsteady techniques 
are virtually identical. In order to solve the unsteady problem, the flow is calculated at 
successive intervals of time by exercising the steady procedure step by s tep m ^time st arting 
from some initial airfoil position and proceeding along the airfoil flight ^patff Th repre 
sentation of the wake can be as simple as a senes of point vortices spread along a straight 
line or as sophisticated as continuously distributed vorticity along a curved sheet, h 

of the convection of fluid particles being attached to it. At each time 
interval the airfoil supplies the wake with some more vorticity (provided the circulation 
around the airfoil is changing), and the already existing vorticity ^ moves ^ rresp °|J^^ - 
downstream. Only this mechanism is essential for the successful outcome of any techniq, 
while the choice of singularities, the use of higher order panels and singulanty distributions, 
as well as the sophistication of the wake model are of secondary importance. 

Regarding the numerical approach, in particular how the airfoil is 

and which tvpe of singularity distributions is used, the present approach follows closely he 
original panel method of Smith and Hess /24/, while with regard to the modelmg of the 
wake, the present approach adopts the procedure as advocated by Basu and Hancock /*/. 
In detail- The airfoil surface is divided into jV planar elements, termed panels, such that the 
airfoil contour is approximated by an inscribed polygon. The first element (index /) starts 
on the lower surfacT at the trailing edge, then elements proceed clockwise ^ound^e airfod 
contour, until the last element (index iV) reaches again the trailing edge .A unifo™ source 
distribution and a uniform vorticity distribution is placed on each panel whereby the source 
strength {a.) k varies from element to element, while the vortex strength y* is the same for 
all elements at a eiven instant of time. The wake consists of a single vorticity panel attached 
as an additional element to the trailing edge, through which the vorticity is shed into the 
wake, and a series of discrete point vortices which are earned dow-nstream with the particles 
of the surrounding fluid. A uniform vorticity distribution of strength (y J* is placed on the 
wake panel, which is further characterized by its length A* and its inclination 0* wi h re- 
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spect to the x-axis of the airfoil frame. After each time step the vorticitv of the wake panel 
is assumed to be concentrated into a single point vortex, which detaches from the airfoil 
and starts moving with the stream. At the same time one wake element detaches from the 
airfoil, a new one is created. The downstream wake of point vortices is thus formed bv the 
shed vorticitv of previous timesteps. 

The following considerations assume that the calculation of (k-I) timesteps has been 
completed and refer consequently to the solution of the kih timestep. At each timestep the 
location of all singularities is frozen and hence any relative motion between singularities is 
not accounted for. The boundary condition of tangential flow is imposed at collocation 
points, which are taken as the exterior midpoints of the elements. Using the concept of in- 
fluence coefficients the tangential flow condition of the ith element can be written according 

jl A ij( a j)k + VkTjBy + (Ks+\)k (Yw)k + Z (Qm)i (r m _! - T m ) = 

7=1 JmX "<=' _ ( 2 . 7 ) 

[v A + 

This condition equilibrates the normal velocity of a fluid particle, that an observer at the 
midpoint of the ah element would notice in undisturbed flow (right hand side) and the 
normal velocity at the very same place due to the perturbation of the flow (left hand side). 
Ihe following singularities are assumed to contribute to the flow perturbation: the source 
(first tom) and vorticitv distributions (second term) on the airfoil contour, the vorticitv 
distribution of the wake panel (third term), and the discrete point vortices in the wake 
(fourth term). The coefficients A, B, and C represent the effect from one panel (point 
vortex) to another panel, because of which they are called influence coefficients. More 
precisely, they are defined as the velocity induced by the singularitv distribution of a single 
element (or by a single point vortex) of unit strength. The first subscript of an influence 
coefficient identifies here the location where the velocity is induced, the second subscript 
denotes the element (or the point vortex) which is causing the very velocity, and the 
superscript indicates which component of the induced velocity is given by the denoted co- 
efficient. The above introduced influence coefficients mean thus the following: 

• A" } represents the normal velocity component at the midpoint of the ith element due to 
a source distribution of unit strength on the jth element. 

• 5" represents the normal velocity component at the midpoint of the ith element due to 
a vorticitv distribution of unit strength on the jth element. 

• Qm represents the normal velocity component at the midpoint of the ith element due 
to the mth point vortex of unit strength. 

The auxiliary boundary condition, which is supposed to uniquely establish the flow sol- 
ution, is taken as the Kutta condition of zero load at the trailing edge (p = p, ), The 
kutta condition is satisfied in an approximate fashion in that the pressures of upper and 
ower surface are equated at the midpoints of the elements adjacent to the trailing edge 
(first and last panel). Relating the rate of change of the difference in the velocitv potentials 
of upper and lower surface at the trailing edge with the rate of change in the* circulation 
around the airfoil, the Kutta condition can be expressed in terms of the vorticitv strengths 
of the current and the previous timestep ' * c ” 


[(» 5)*] 2 - Ov )*] 2 = 2 [ 


f d{4>n~4>\) 1 

- 2 r -^-i 

L dt J 

k 2 L dt J, 
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where T denotes the circulation around the airfoil and / is the airfoil perimeter. The time- 
derivative of the circulation is here approximated by a first-order accurate backward finite 
difference. 

In contrast to the steady flow solution there is one additional condition to be imposed on 
unsteady flows. The Helmholtz theorem postulates that any change in the circulation 
around the airfoil must be matched by a change of vorticity in the wake of equal magnitude 
and opposite sign. If this postulate is evaluated in terms of our wake model, the following 
relationship is obtained 

A* (yj k = l (y k -\ ~ Vfc) ( 2 - 9 ) 

This equation states that the vorticity in the wake element be equal to the negative change 
of circulation around the airfoil with respect to the previous timeline. This equation is fur- 
ther based on the assumption that the circulation around the free point vortices in the 
dow'nstream wake does not change after the point vortex is formed from the spread 
vorticity along the wake panel. 

Before we continue with the actual algorithm, we summarize briefly the setup. The follow- 
ing unknowns have been introduced in the equations (2.7) through (2.9): There are the <V 
unknown source strengths (<r ; )* and the one unknown vorticity strength y k along the ele- 
ments on the airfoil, and there are the unknown vorticity strength (yj*, length A* , and 
orientation 0* of the wake panel. The last two of which did not show up in the equations, 
but they have to be known as part of the geometrical setup in order to calculate the influ- 
ence coefficients due to the wake panel. This makes a total of A -I- 4 unknowns, which must 
be matched by an equal number of equations. So far the following equations have been 
identified: 

1. The flow tangency conditions represent a system of N equations. 

2. The Kutta condition is a single equation that determines the circulation around the 

airfoil. . , 

3. The Helmholtz theorem is another single equation that provides a relation for the 
strength of the vorticity distribution along the wake element. 

Clearly we are short of two equations, which ought to address the length and orientation 
of the wake element. Following a suggestion by Basu and Hancock, the following assump- 
tions are made about the geometry of the wake panel: 

1. The wake panel is oriented in the direction of the local resultant velocity at its midpoint 

tan ©»--)% (2-10) 

with (o*) A and (t^)* denoting the x- and y-velocity components at the midpoint of the 
wake panel, as viewed in the (moving) airfoil frame of reference. 

2. The length of the wake panel is proportional to the magnitude of the local resultant 
velocity at its midpoint and to the timestep 

A* = yj [(^w)*] 2 + [(°w)itP ( [ k~ l k- 1) (2-11) 

w’hereas the velocity components (u*)* and (i/J k (in both equations (2.10) and (2.11)) 
do not include the effect of the wake element on itself. 

Although the majority of the equations is linear, the nonlinearities in the Kutta condition 
and in the tw r o equations for A k and Q k necessitate an iterative solution procedure. How- 
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ever, the effect of the nonlinearly involved variables (A* and 0*) is weak enough to allow 
for a linear algorithm being applied to the system of flow tangency conditions, in which the 
nonlinear variables are treated as known quantities given by the previous cycle of iteration. 
Each cycle closes by recomputing improved values for the nonlinear variables. Thus the 
overall solution is obtained in successive approximations, which involve the following steps: 

1. Before the calculation of a new timeline can get started, the length and the orientation 
of the wake element must be guessed, for example from the previous timeline. 

2. Making use of the Helmholtz theorem (the vorticity strength of the wake element is 
replaced by an expression of the unknown y k and the known y A _,) the system of flow 
tangency conditions can be rearranged in form of the following matrix equation 

[^] {°)k ~ Vk{B)k + {Q* (2.12) 

where [ A ] is a N x N matrix, {5}* and {C} 4 are N-dimensional vectors, respectively. 
This matrix equation represents a linear system for the unknown source strengths (crd* 
in terms of the yet undetermined vorticity strength y k . 

3. With these results the Kutta condition turns into a quadratic equation for the unknown 
vorticity strength y k . 

4. Under the assumption of an accurately determined wake element steps 2 and 3 establish 
a unique solution of the flow. Now the length and the orientation of the wake panel 
are recomputed and the procedure proceeds with step 2 if A k and Q k have not con- 
verged to a desired accuracy. 

5. In case of convergence the calculation of the kth timeline is completed by determining 
the velocity and pressure distributions on the airfoil. 


Due to the minor effect of an inaccurately assumed wake panel on the overall flow results, 
the procedure converges extremly fast, say within two and four iterations. The numerical 
procedure can further be accelerated by taking account of the fact that the matrix [/l] does 
not change with the iterations nor with time. Hence the most time consuming numerical 
step, which is the LU -decomposition of the matrix [/f], needs to be exercised onlv once. 
However, a record of the manipulations performed on the right hand side during the initial 
timestep must be kept in order to repeat these manipulations on the right hand sides of the 
following iterations and timelines. 


Once the source strengths and circulation around the airfoil have been determined, one can 
proceed with the computation of the velocities throughout the flowfield. Utilizing again the 
concept of influence coefficients the tangential velocities at the midpoint of the ith element 
can be expressed in the following form 

n ,v k — i 

(*?)* = 1.4 (»;)* + Yk ZsJ + (BU +l )k(vJk + Z(cU(r m _ t - rj (2.13) 

V— i y=l m=\ 




(2.14) 


{Vf) k in equation (2.13) refers here to the tangential velocity as viewed in the inertial frame 
of reference, whereas (of)* in equation (2.14) refers to the tangential velocity as viewed in 
the moving frame of reference. If there is interest in the load on the airfoil, then the un- 
steady Bernoulli equation must be resolved with respect to the pressure to provide this kind 
of information. The pressure coefficient in an unsteady flow is given by 


<=>= - 


2 

tj 2 dt 


+ 


• r2 


(2.15) 


10 


2. The unsteady potential flow method 


In contrast to the steady flow, the unsteady pressure coefficient includes a contribution due 
to the rate of change of the velocity potential. Since the value of .3 cp/dt will be approximated 
bv a first-order accurate backward finite difference, the velocity potential must be evaluated 
at each timeline. The calculation of the velocity potential is accomplished by an integration 
of the velocity flowfield along a line extending from a position, which is sufficiently far 
upstream so that disturbance effects there are negligibly small, to the point of interest. In 
particular, the integration is performed in two steps, the first of which proceeds along a 
straight line from some upstream position to the leading edge of the airfoil, and the latter 
continues then along the airfoil contour from the leading edge to the midpoint of each 
panel. 

The setup for the next timeline requires the convection of the singularities within the wake. 
The free point vortices of the wake are carried with the fluid as if they were fluid particles. 
The distributed vorticity along the wake element is first concentrated into a single point 
vortex, and then converted from the midpoint of the wake element in the same lashion. I he 
path of each particle is obtained by a single step integration procedure (predictor step only), 
which approximates the actual path by short pieces of straight lines that are tangent to the 

actual path 

jX.( r fc+i)j = + ( V m ) k ( t k+l ~t k ) ( 2 - 16 ) 

This completes the formulation and description of the method to solve the in viscid part : of 
the flow over an airfoil that is executing some kind of time-dependent motion, for further 
details, like the calculation of the influence coefficients, the reader is referred to Teng (thesis 
at NPS 1421). 
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3. The unsteady boundary layer method 

Introduction 

For many flows of aerodynamic interest the Reynolds number is large enough such that the 
flowfield can be divided into two regions: an inner region, which consists of the boundary- 
layers along the wall and in the wake, and an outer in viscid region. The flow within the 
inner region is assumed to be governed by Prandtl s boundary layer equations, which con- 
stitute one of the asymptotic approximations to the Navier-Stokes equations for large 
Rey nolds numbers. The Navier-Stokes equations are thought of as the exact equations of 
motion for a viscous fluid. Furthermore the full equations (which are written in instanta- 
neous flow quantities) apply equally to laminar and turbulent flows. Though the full 
Navier-Stokes equations would resolve every detail of the turbulent motion (no matter how 
small the scale is chosen), their solution is, at the present time, not seriously considered for 
engineering purposes because of costs. In order to solve turbulent flows for applied prob- 
lems the time-accurate equations are replaced by time-averaged equations. This is accom- 
plished by time-averaging the turbulent fluctuations over a time-scale that does not resolve 
any details of the turbulent motion, yet is small enough to describe every 7 detail of the mean 
flowfield. Virtually all engineering methods pursue this approach, that is the solution of the 
Reynolds-averaged Navier-Stokes equations or of some asymptotic approximation to them. 
It is owing to Prandtl, that the Navier-Stokes equations had been reduced to a much sim- 
pler form opening the way to a wide range of methods, that enable to account for viscous 
effects in a simple fashion. Prandtl's boundary layer concept derives from the assumption 
that changes in the velocity (within the viscous layer) occur more rapidly across the layer 
than along it. The larger gradient of the velocity in the cross-stream direction can be ex- 
plained by means of the no-slip condition, which requires that fluid particles on the airfoil 
contour must not move with respect to it. Consequently, there is only a thin layer left, in 
which the particles shall change their velocity from zero to the external frictionless velocity. 
Such a layer is referred to as a boundary layer, or (using the more recent terminologv) as 
a thin shear layer. 

The nature of the flow in a boundary layer differs considerably from that of Navier-Stokes 
governed flows. Whereas the Navier-Stokes equations provide a mechanism such that dis- 
turbances get carried in all directions, the boundary layer equations put some limitations 
on that process. Foremost, the streamwise diffusion, asencountered in Navier-Stokes gov- 
erned flows due to viscous and turbulent stress gradients, does not take place in boundary- 
layers. Furthermore, there are no pressure signals propagating within the boundary layer, 
because the pressure is imposed from outside. This gives rise to a flowfield, in which the 
flow is uninfluenced by downstream events. Such a behavior is, mathematically speaking, 
termed parabolic. From a numerical point of view, the parabolic behavior is one of the 
biggest assets of the boundary layer equations, since parabolic equations admit marching 
algorithms. These algorithms obtain a solution by proceeding step by step in the marching 
direction such that the solution at a particular point is made dependent on the results of 
already computed steps only'. The unsteady boundary layer equations offer two of these 
directions, namely the streamwise direction and time. The presented methodologv takes 
advantage of both: For a given time level the procedure performs a sweep for each the up- 
per and the lower surface, with the solution of each sweep advancing step by step from the 
stagnation point towards the trailing edge. After a complete solution of a particular time 
level has been obtained, the scheme proceeds to the next time level and the described rou- 
tine is repeated. Accordingly, a single step will yield the solution for a column of points at 
a given streamwise position and for a given time. In contrast to the stepwise advancing of 
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the solution in the streamwise direction and in time, the solution for such a column must 
be attempted simultaneously in order to account for the viscous diffusion in the cross- 
stream direction. Each step adds another column to the already computed flow solutions 
such that the field (x, t-plane) is successively covered by those columns. The efficiency of 
the boundary layer approach derives from this algorithm of marching solutions, which by- 
passes the laborious field solution as required by the Navier-Stokes equations. 

While the boundary layer approach represents a very economic means of computing 
viscous flows, it is not a universal tool. The boundary layer equations do apply only to 
flowfields, which satisfy the assumption on which the boundary layer approximation is 
based. That is the condition that viscosity affects the flow essentially only in a very thin 
laver. This assumption excludes flowfields experiencing a massive separation, like the flows 
past blunt bodies or the flows past streamlined bodies at high angles of incidence. Further 
because the pressure is imposed on the boundary layer by the outer flow, the classical 
boundary' layer approach is not applicable in regions, where the viscous layer has a sub- 
stantial effect on the pressure. As far as the boundary layer is concerned, the pressure is a 
known function of the streamwise coordinate and of time. The pressure ceases to be a 
function of the cross-stream coordinate and hence remains constant across the boundary- 
layer. Because of this the boundary layer approach cannot cope with flows past bodies, 
whose surface curvature either is large or changes abruptly. In spite of all these limitations 
Prandtl's boundary' laver concept proves extremely useful in analysing viscous flows be- 
cause of the substantial simplifications achieved by reducing the originally governing 
Navier-Stokes equations. 


The differential equations, the boundary and initial conditions 

The boundary layer equations can be obtained by applying an order of magnitude analysis 
to the Navier-Stokes equations. Utilizing Prandtl s finding, that the thickness of the viscous 
layer is small compared to a characteristic length in the primary flow direction, the magni- 
tudes of the terms in the Navier-Stokes equations can be estimated. Keeping only those 
terms of the largest order of magnitude leads to Prandtl s famous set of boundary la\er 
equations (see for example / 37/). The time-dependent boundary layer equations for two- 
dimensional, incompressible flow take in terms of dimensional variables the following form 

^L + 4^ = 0 (3-1) 

dx dy 


cu 

dt 


+ u 


ou 

dx 


+ o 


cu 

dy 
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8u e 
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oy 




(3.2) 


The foregoing equations are subject to the usual boundary conditions, which are the no-slip 
condition on the airfoil contour and the prescription of the external velocity at the outer 
edge of the boundary layer 


j; = 0: « = D = 0 

y — * oot u * u e (jc, /) 


(3.3) 


Here the b denotes a scaled viscosity such that the equations apply equally to the laminar 
and turbulent regimes. Its turbulence related contribution must be provided by a turbulence 
model, which, in our scheme, is based upon the eddy viscosity concept. All flow variables 
represent ensemble averages over a large enough time, such that the above equations need 
not describe the turbulent fluctuations. The equations are expressed in curvilinear coordi- 
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nates x and y, which are directed along and normal to the airfoil contour, and in the ve- 
locity components u and o along these directions. The pressure gradient in .r-direction is 
written in terms of the tangential velocity component of the external flow, which, to a first 
approximation, may be assumed identical to the frictionless flow past the airfoil. Conse- 
quently the external velocity distribution is obtained from an unsteady panel method, which 
solves for the (inviscid) potential flow past the airfoil. 

In order to complete the formulation of the problem, initial conditions must be prescribed 
for the entire .r.y-plane at some initial time as well as upstream boundary conditions for 
the entire f,y-plane at some initial section x s 

r=t t and x> x s : u = u i (x,y) 

x = x s and t>t i \ u = u s (y,t) 

On the assumption that steady flow conditions prevail at the time the initial conditions 
can be obtained from a steady-state solution of the boundary layer equations. Further, it 
is assumed that the unsteady motion develops gradually (dujdt = 0 at r,), such that double 
structured boundary layers need not be considered. The upstream boundary conditions re- 
quire the prescription of a velocity profile at the stagnation point for r > t.. While the gen- 
eration of the upstream boundary conditions does not present any problems "for 
stead} -state flows (recall that the velocity profile at the steady stagnation point can be de- 
rived from a similarity solution), this undertaking proves much more difficult in the case 
of unsteady flows. First, because the tangential velocities at the stagnation point in un- 
steady flows do no longer vanish across the layer, as they do in steady flows, and second, 
because of the possible occurrence of flow reversals at the unsteady stagnation point. The 
details of the procedure that can cope with both problems are given in Appendix C. In 
summary, the procedure consists of an iterative computation of the flow at the stagnation 
point, and the adjacent upper and lower surface points. In particular, the computation of 
the flow at the stagnation point involves the use of a modified characteristic box scheme, 
which differs in that the momentum and continuity equations are solved successively and 
in a decoupled manner. In closing this paragraph it seems appropriate to point out", that 
only differential equations' that are accompanied by suitable initial and boundarv conditions 
constitute a mathematically well posed problem. 


Signal propagation in unsteady boundary layers 

The discretization of the just established, unsteady boundary layer problem shall be de- 
scribed next. Formally, the partial differential equations are approximated by finite differ- 
ence equations, w-hich can be found through the use of Taylor-series expansions. Though 
such an approximation will be consistent with the exact equations (in the sense that in the 
limit of an infinitely small meshsize the truncation error of the approximation tends towards 
zero), it just might not work. The necessary - and in case of a linear equation sufficient - 
condition for a numerically converging scheme is stability. Stability means here the capa- 
bility of the scheme to dampen errors of any source while advancing from one marching 
step to the next. From a physical point of view even such a scheme might be faultv. A 
physically sound scheme needs, in addition, to incorporate the mechanism of signal propa- 
gation as imposed by the differential equation. In the following we will elaborate on this 
process and its consequences on the design of a finite difference'method. 

The unsteady boundary layer equation is characterized by two well known tvpes of signal 
propagation: Regarding the streamwise velocity u as a function of x,y as well as ofy, t the 
equation is diffusive, whereas regarding u as a function of x , i it behaves like a wave 
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equation. Consequently, the unsteady boundary layer equation contains two time-like 
initial-value variables, 'namely the streamwise coordinate * and the time /, and one 
boundarv-value space variable, namely the boundary layer coordinate y. How does this re- 
late to the mechanism through which information is transmitted in an unsteady boundary 
laver? Let us assume that a slight disturbance is generated at the point (x lt y { , fl) and let us 
trace its route. In the very instant the disturbance is generated, it propagates immediately 
across the boundary layer station (jq, /,), that is along the line * = jq and / = /, all the way 
down to the wall and up to the edge of the boundary layer. In a local perspective lnlorma- 
tion also travels along the streamlines with the speed of convection. But globally this 
process adds up to a streamwise transport of information with max(u) and iran(u) as the 
velocities of propagation. This is readily explained. The disturbance instantly diffuses to 
all lavers of the station (*„ /,), from where it continues its path in each of these lasers 
0 = constant planes) with the local particle speed u. As soon as the signal reaches a new 
station (x-,,y 2 , z 2 ), the signal there immediately spreads to all values of y at that station. 
Because of this mechanism the signal effectively travels with the maximum velocity m the 
downstream direction, and in case of a Bow reversal with the maximum backflow velocity 
in the upstream direction. Thus the mode of propagation in the streamwise direction is as 
follows: The signal travels up to the edge of the boundary layer in the instant of its gen- 
eration, and then travels through the outmost part of the boundary layer with the speed 
of the external flow, which carries the signal downstream at the fastest speed possible. At 
each new station the signal is felt instantaneously across the layer because of the viscous 
diffusion in y-direction. There are no signals propagating upstream if there are no flow re- 
versals. In case there are, signals follow a pattern similar to the downstream propagation, 
with the maximum reversed flow velocity taking the role of the external flow velocity. 

As a consequence of the above outlined mechanism, the flow at a particular point does not 
depend on the entire flowfleld, not even on the entire upstream llowfield. The formal tool 
to describe the flow of information is the concept of the zones of dependence and influence. 
They are defined as the zone whose points - and only w'hose points - do affect the now at 
a particular point, and as the zone whose points - and only whose points - ca ^ be reached 
by signals sent from a particular point. For the unsteady boundary' layer equations th ^ z ° n 
of dependence (zone of influence) is a wedge-shaped volume (in the x,y, /-space), which is 
generated bv two planes that are perpendicular to the x, /-plane and that include the out- 
ermost directions of all particle paths (projected in x, /-planes) across the boimdary lavei. 
These two directions correspond to the maximum streamwise velocity, which is usually 
external velocity, and the minimum streamwise velocity, W'hich is either zero or the largest 
reversed flow' velocity 


= min [u(x,y, /)] = — max |u rev l 
dt 0<y<<5 0 <y<6 

= max [u(x,y, /)] = maxjul 

dt 0<y<<5 0<y<<5 


(= 0 ) 
( = ««) 


(3.5) 


The implications of these findings regarding the layout of the finite difference scheme 
usually summarized in two rules: the dependence rule and the law of forbidden signals. The 
former one requires that the computational domain of dependence brackets the domain o 
dependence as defined by the differential equation. This rule is occasionally linked to the 
ubiquitous Courant- Friedrichs- Lewy condition / 15/, which as it is argued, can be inter- 
preted geometrically in such a way. While that is true for the original CFL-condition, which 
is the stability criterion of explicit schemes for solving linear hyperbolic equations, the 
generalization to parabolic problems and implicit schemes does not hold good. In general, 
the stability constraint may not coincide with the condition resulting from the dependence 
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rule. In fact, schemes have been reported that do satisfy the dependence rule and are in- 
herently unstable /2S/. The second rule, the law of forbidden signals, is rarely recognized 
and it is violated m most schemes currently used. This law postulates that numerical 
sc he m es shall prohibit the reception of signals that cannot exercise anv influence, because 
tv? 7 d ,u n ? onginate W1 , thjn the zone of dependence as defined by the differential equation. 
™ us ! h< r de P endence .rule and the law of forbidden signals lead to contradictory statements 
though both are motivated by the proper implementation of the physics. Well there is a 
scheme possible that can satisfy both rules, namely one in which the computational domain 
of dependence and that one defined by the differential equation do coincide. For now the 
objective is to find a scheme that implements these principles as well as possible under the 
constraint of stability and in consideration of the limited computer resources available. 


The numerical procedure 

The numerical solution of the unsteady boundary layer equations is obtained with a tech- 
mque that utilizes both the regular box scheme and the characteristic box scheme. The box 
method was devised by Keller /2 5/ for solving diffusion problems, and has, since then, ex- 

PP ff d t0 3 m de Va A 1Cty °r P arabolic Problems /26/, with particular success 
/ S 5 ? n ° w Problems. One of the strongest advocates of the box method has 

S®!? ls) T C . ebeC1 ’ Who ’ together with his collaborators /3,10,1 1/, has contributed a stood 
de a i of the techniques associated with this method. The differencing of the box method is 
implicit with respect to the cross-stream direction, and second order accuracy is formally 
^ranted m time and both space coordinates on arbitrary, nonuniform nets. The' box method 
is said to be unconditionally stable, even with the use of a large meshsize and the appear- 
ance of rapid net variations. w 

After this brief introduction of the box method, we focus on the peculiarities of the un- 
steady-approach, in particular in comparison with the well-established steady box method 
1 hough the unsteady boundary' layer equations differ from their steady counterparts onlv 
by the appearance of the time-derivative du/dt, the concepts of solution for steadv flows 
o not directly carry over to unsteady flows. This is not so much because of the additional 
dimension, rather than because of the different nature of the equations, in this regard the 
unsteady, 2-dimensional equations are much closer related to the steady, 3-dimensional 
equations The thread between the two of these is the common mechanism of signal prop- 
agation, which, in turn, is instrumental to the choice of the finite difference representation. 

e essential feature, missing in practically all steady, 2-dimensional procedures, is the for 
unsteady, ^-dimensional and steady, 3-dimensional flows absolutely crucial, incorporation 
of upstream influence. The presented technique takes this fully into account bv switching 
‘T * he 5 eg ? lar t0 .l be characteristic box scheme /8/ as soon as the streamwise velocity 
indicates backflow;. The characteristic box scheme admits the indispensable upstream in- 
fluence by expressing the governing equations with respect to those directions, along which 
are df 311 ?^ 111 ^ /_P ane are locall >' Propagated. The integral curves of these directions 


dx(s) di{s) 

~7T = u ~dT ( 3 - 6 ) 


where s denotes the coordinate along such a curve. The method owes its name to the con- 
cept of characteristics, which is well knowm in the theory of hyperbolic equations but also 

a S 1nnau-h te L y i mke t t0 the m< : cha ™ sm of si S nal propagation. The abov e q addressed paths, 
along v Tuch disturbances are locally convected, are associated wdth so-called subcharacter- 
istics. (These are characteristics of a simplified problem, in which the highest derivatives of 
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the original problem, in our case the viscous diffusion term, are neglected. For a more 
elaborate explanation see for example reference /40/.) The characteristic method uses then 
differencing along the subcharacteristic direction 

5 + M 8 _ _d_ 

dt dx dt 


s^const 


1 d 

dtisMds ds 


(3.7) 


such that the momentum equation can be reformulated as 


1 dlL = _ 0 l!L + l^ + u Jh. + v J-(b^-) (3.8) 

dt(s)jds ds cy dt dx ay v dy 

The significance of this equation is that its solution corresponds formally to an integration 
of the "right hand side along the subcharacteristic curve. This circumvents the difficulties 
encountered when the regular method is employed in regions of reversed flow. 

Before we continue with the details of the box method, we introduce a coordinate trans- 
formation. Most of the boundary layer methods employ scaled variables in order to reduce 
the rapid variation that occurs near the wall. Numerically it is advantageous to deal with 
quantities of the same order of magnitude. To this end the dependent variables are 
nondimensionalized by a reference length L and a reference speed U x , which in the actual 
calculation are taken as the chordlength and the freestream velocity, respectively 


i = 





(3.9) 


where the Reynolds number R L is defined in terms of the reference scales 

R L = -^~ ( 3 . 10 ) 

V 

Not accidentally the normal coordinate is stretched by a factor involving the squareroot 
of the Reynolds" number. It is one of the first and simplest_results of boundary layer theory, 
that the boundary layer thickness is proportional to 1 Id Rl ■ Thus the newly defined, normal 
coordinate y will’ be "of the same order as the streamwise coordinate. This particular coor- 
dinate transformation (yielding so-called physically transformed coordinates) is not as 
powerful as other commonly used transformations, say for example the Falkner-Skan 
transformation. It does neither maintain a constant boundary ;ayer thickness nor eliminate 
the stagnation point singularity, as the Falkner-Skan transformation is able to. The reason 
for introducing such a simple transformation lies to some extent in the necessities dictated 
by the unsteady formulation of the box method. 

The solution of the boundary' layer equations by the box method involves four steps, of 
which the first is discussed later in this section, the second and third are described in Ap- 
pendices A and B, and the fourth must be looked up elsewhere ( for example in any textbook 
by Cebeci / 3,6/). The four steps are: 

1. The boundary layer equations are reduced to a system of first order differential 

equations. . 

2. The first order differential equations are approximated by simple centered differences 
and two-point averages, using values at the comers of one difference molecule only. 

3. The resulting algebraic equations are linearized by Newton's method. 

4. Finally, the numerical solution of the resulting block tndiagonal system is carried out 
by the block elimination method. 
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The first step of the box method is in order to allow the approximation of all functions and 
derivatives in terms of the quantities at the comers of one difference molecule only. To this 
end a new dependent variable, V, is introduced to eliminate the second order derivatives in 
the momentum equation. Using the scaled velocities U and F, the scaled vorticity V, as well 
as the scaled external velocity W 


U = 


U 


v= 


L 8 u 


JR l ^oo dy 
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W = 




(3.11) 


the boundary layer equations can be rewritten as a system of three first order differential 
equations 
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(3.12) 
(3. 13. a) 
(3.13.6) 
(3.14) 


The first of these equations refers to the conservation of mass, the second and the third 
represent the conservation of momentum either in regular or in characteristic form, and the 
fourth defines the extra variable. Equations (3. 13. a) and (3.13.6) are used alternately, de- 
pending on the sign of the local streamwise velocity. The foregoing equations are comple- 
mented by boundary conditions, requiring no slip at the wall and matching of a prescribed 
velocity at the edge of the boundary layer 


>7 = 0: U = F = 0 

V ~ Ve : U=W{$, r) 


(3.15) 


When completed by a closure relationship for the turbulent shear stresses (see following 
section), by appropriate initial and upstream boundary conditions (see preceding section 
and Appendix C), the problem is all set for the numerical solution. For the details of the 
further proceeding the reader is referred to Appendix A for the regular and Appendix B for 
the characteristic box method. 


Before we close this section, we would like to comment on a practical and a hypothetical 
aspect of the solution procedure. The practical implementation of the regular box scheme 
required a modification of the commonly used procedure. The modification was necessi- 
tated by the occurrence of streamwise oscillations in the near-wall solutions. An analysis 
of this instability, together with the outline of the modification, is given in Appendix E. The 
comment on the hypothetical aspect deals with the question how well the box method sat- 
isfies the dependence rule and the law of forbidden signals. Following the explanations of 
the preceding section the domain of dependence is a wedge-shaped volume, which is gen- 
erated by two planes that touch the outermost subcharacteristics. The dependence rule re- 
quires that information be included from all points located between those two planes. While 
the characteristic box method attempts to simulate the local mechanism of signal propa- 
gation, and it succeeds in that by formulating the momentum equation along'the path a 
disturbance would take in a local perspective, it fails to reproduce the global mechanism. 
The characteristic box scheme includes information from a single point only, namely the 
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origin of the subcharacteristic on the previous timeline, whereas all points between the two 
outermost subcharacteristics should be included. Consequently the characteristic box 
method is seen to violate the dependence rule, though not to a dramatic extent. On the 
other hand, the characteristic scheme certainly satisfies the law of forbidden signals. Nota- 
bly the regular box method does not. Information on the same timeline, that is included 
in the regular scheme, is not only undue, but physically impossible since signals do not 
spread with infinite speed in the streamwise direction. Finally, it should be pointed out, that 
the technique of regular and characteristic boxes probably reproduces the mechanism o 
signal propagation better than any other method. 


Turbulence modeling and transition 

In order to solve the boundary layer equations for the turbulent flow regime, the governing 
differential equations must be completed by a closure relationship for the apparent tur- 
lent shear stress. This term has been added to the equations in the process of averaging 
time-accurate equations, which is in order to bypass the extremely costly calculation of the 
turbulent fluctuations. Furthermore, from the engineers point of view the compu a 1 
of the turbulent motion beyond a certain time- and length-scale seems unnecessary, since 
the main focus is on the mean flow. The averaged equations combined with a turbulence 
model represent then, hopefully, a less complicated set of equations. Turbulence modeling 

consists of providing suitable means for the approximation of the turbule ! lt t h e eddv 
Most of todav's models fall into two categories: those which are based upon the eddv 
viscosity hypothesis, and those which solve one or more Reynolds stress transport 
equations. Both tvpes of models involve empirical constants and hence require a venfica- 
tiTn bv comparison with experimental data. Establishing a relationship between the un- 
known turbulent shear stresses and the mean flow variables closes the set 
Reynolds-averaged equations in that the number of unknowns is reduced to equal the 

number of equations. 

Todav there is a wide variety of turbulence models available, ranging from the simple al- 
gebraic eddv viscosity models to the sophisticated stress-equation models. The classifica ion 
of models is according to the number of the partial differential equations . that .are ^olved 
in the calculation of the shear stresses. Though the more complex models, including 
widely used k-i model, promise more generality, it is not obvious that their accuracv is su- 
perior to simple models. The most common approach (in engineering calculations) adopts 
Boussinesq's concept of eddy viscosity ( 1877), in which the apparent turbulent shear stress 
is assumed to be related to the rate of mean strain through a turbulent (eddv) viscositv 


r 

P 


X, + T, 


= V 


du 

dy 


— u'v' = (v -t- v f) = vbjfr 


(3.16) 


This model simplv adds a turbulence related viscosity to the laminar viscosity. T e 
strengths of the concept are its simplicity and the convenience of maintaming the same 
form of differential equations for laminar and turbulent flows. The major shortcomings of 
these models lie in that the turbulent viscosity is evaluated in terms of local flow parameters 
and in the need to separately tune the model for different classes of flows. In view of the 
engineering oriented purpose of our efforts, we decided m favor of such asi mpemo d el.m 
particular we use the algebraic eddy viscosity formulation of Cebeci and s ^h /6/. This 
model, representing a zero-equation closure, evaluates the turbulent shear stress through 
an algebraic equation. The eddy viscosity is given by a two-layer model, whose inner region 
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is based on Van Driest s approach, and whose outer region is based on a velocity defect 
approach 


V ' - {o.4j-[l -exp(-i )]} J | M 
V' = 0.0168 If {tig u) dy | y (r 


for 0 < j; < ^ 
for y c <,y<,8 


(3.17) 


where A = 26v/(v )' /2 

' Oy 'max 

Here the crossover distance^ is defined as the location where the eddy viscosities of inner 
and outer regions coincide. The presented model is essentially identical to the steadv-state 
formulation. Originally, Cebeci and Smith did include some unsteady correction terms, 
which showed up in the damping constant A. In a previous investigation, however, the 
unsteady correction terms turned out to have little effect on the results /20/. Consequently 
in our computations these terms are not accounted for. Nevertheless, there is plenty of ef- 
fects that require alterations to algebraic turbulence models. For example it is weirknown 
that transition is not an instantaneous process, and hence the transitional region ought not 
to be reduced to a switching point between the laminar and the turbulent regions. In order 
to account for the transitional flow region, we make use of the intermittencv distribution 
by Chen and Thyson /14/ 


Vir - 1 - exp 






(3.18) 


where x„. and are the location of the onset of transition and the transition Reynolds 
number, respectively. The empirical factor G y exercises a major influence on the length 
of transition. The originally assigned value of 1200 is suited well for high Reynolds number 
flows, whereas flows at lower Reynolds numbers require reduced values of G . Modifica- 
tions due to other effects, like surface roughness, freestream turbulence, or streamwise wall 
curvature, have, as of now, not been implemented in our code. Though algebraic eddy- 
viscosity models have accumulated a remarkable record, extreme caution is advised when 

models are applied outside the range of flows, upon whose data the empirical constants are 
based. 


Since flows past airfoils comprise laminar, transitional, and turbulent regimes, there is the 
need to predict the onset of transition. In spite of the unargued importance of the accurate 
prediction of transition, our understanding of that process is unsatisfactory and incomplete. 
The currently available methods include empirical data correlations, which relate the onset 
of transition to some mean flow parameter, and linear stability theory, in which transition 
is deduced from the amplification ratio of disturbances obtained from the solution of the 
Orr-Sommerfeld equation. The prediction of transition in our scheme is based on a combi- 
nation of Michel's empirical data correlation /33/ and Smith and Gamberoni's ^-method. 
According to Cebeci / 6/ the resulting expression is given by 


\ = 1 - 174(1 
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where R g and R x are the Reynolds numbers based on the momentum thickness and the 
streamwise coordinate at the onset of transition, respectively. 


Unsteady separation and singularities 


The subject of flow separation and the possible occurrence of a singularity incident thereto 
has had a long and controversial history. Separation is the flow phenomenon, which, ac- 
cording to Brown and Stewartson /4/, can be described verbally as follows. Here the main 
stream 8 which has hitherto been in close contact with the body, suddenly, and for no obvi- 
ous reason breaks awa-y, and downstream a region of eddying flow, is set up. i 
Brown and’stewartson continue to quantify their description in that they define separation 
as the transition from an attached boundary layer which confines the viscous effects _to a 
laver of 0(R~ li2 ), to a detached boundary layer, where viscous effects are felt in la ers 
ther from the wall. However, neither of these definitions is suited for identifying the point 
of separation in practice. Thus a considerable number of separation criteria has been pro- 
posed, but their general applicability, including time-dependent and 3-dimensional flows, 
remains to be proven. 

The question of ample importance is whether flow separation can be related to the occur- 
rence^ a singularitv in an uninteracted boundary layer description. For a long time it was 
known that the classical boundary layer solution breaks down due to the apP«"“ c ® ° f 
Goldstein singularity /22 / in steady flows. The singularity can ^ some 

terized bv the nonexistence or nonumqueness of the solution, or by the growth of t some 
flow variable beyond all bounds. Clearly, the Goldstein singularity is not P h ysical and the 
Navier-Stokes equations have even been shown to be free of such a singulan ty ■■TJje singu- 
lar behavior of the solution can rather be read as an exaggeration of the actual flow and 
as such might indicate the need for another, more appropriate set of equations. The exist- 
ence of a similar singularity in the classical, say uninteracted unsteady boundary lay er 
solutions has been the subject of numerous investigations with 

dictorv findings. The controversy was settled by Van Dommelen and Shen /43/, vh g 
a numerical proof that a singularity will develop in the laminar boundary layer of an 
impulsivelv started cylinder within a finite time. Since on the one hand the existence of 
singularities in the classical boundary layer description has been proven, and on the other 
hand the actual flow is not known to exhibit a singular behavior the question arises, why 

a certain description does develop a singularity, and how it can be “ ' Juatfons as 
ance of a singularitv is believed to be the consequence of solving incomplete equations, s 
well as of an incorrect choice of the external velocity distribution. The boundary layer 
simplv must not be considered separately from the external inviscid flow. Instead it is 
necessarv to account for the mutual interaction between the boundary layer and the ex- 
temal flow. This led to the development of the interactive boundary aver met ho s, \v i 
were able to remove the Goldstein singularity in steady flows, and at least to delay the ap- 
pearance of Van Dommelen's singularity in unsteady flows. 

Though the interactive methods enabled the computation of (moderately) separated flows 
bv means of the boundary layer approach, they did not contribute to ^ the yet Q ™ s< ^ e 
of identifying the point of separation. While this is not a problem m steady fiows, where the 
point of separation coincides with the point of zero wall shear, it presents a formidable task 
fn unsteadv nows. There, the point of zero wall shear w,ll. in general, be nerther singular 
nor identical with the point of separation. None of the thus far proposed unsteady sepa 
ration criteria turned out to be easily implementable and at the same time universally ap- 
plicable The generalized Goldstein criterion, which suggests the coincidence of separation 
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r^?n?T S la D c Can 0n ^ y 0b l erve s ^ m P toms ’ and lacks the practicality of a definite corre- 
Th ^ i 1RS ; C u ten ° I ? ^ the a< ; ronym MRS stands for Moore / 34/, Sears /38/, and Rott 
/36/) was deduced by analogy with a steady boundary layer along a moving wall. Unsteady 
paration is predicted at the point, at which the shearing stress and the velocity, with re- 
E a of , refer ^ nce moving with the speed of separation, simultaneously vanish, 
^vr / the doubts about the universality, the implementation of such a criterion seems 
next to impossible due to the a prion unknown speed of separation. The probably most 
promising, though as well not very practical, approach is by Wang /47/. His criterion is 
based on the analogy with the well-established separation criterion for steadv 
3-dimensional flows. It requires the evaluation of so-called skin friction lines (these are in- 
tegral curves whose slope equals the skin friction) in the *, r-plane. Wang then proposes 
unsteady separation at points, at which the skm friction lines coalesce, provided thev do so 
at all. Because of the lack of experimental evidence the status of any unsteady separation 
criterion remains for the time being, that of uncertainty. Obviously/further work, both of 

experimental and numerical nature, is needed to attain a clearer picture of unsteady flow 
separation. * 


22 


3. The unsteady boundary layer method 


4. Results and discussion 


Using the technique discussed in the preceding sections, pilot calculations were performed 
in order to do a system checkout prior to developing the interactive version of the code. 
The reader is reminded of the interim status of this investigation, which permits to run the 
code in the direct mode only. For that reason the tedious procedure of a detailled compar- 
ison with experimental data has been saved for the final version of the code The calcu- 
lations conducted thus far refer to the unsteady flow past an airfoil undergoing a 
ramp-change in angle of attack. Results were obtained for a small series of cases different 
in the Revnolds number, the initial and final angle of attack, and the time required for the 
pitch-up/The presented results are confined to a single case, which can be described as 
follows. A NACA 0012 airfoil changes its angle of attack according to the rule 


for r < 0 


oc(t) = 


ft . i'Jl t V 

f 1 

) or, + (ay - a,) (3 - 2t/ T/ ) t 2 /t) for 0 < r ^ y 
l cut f° r T > T / 


(4.1) 


For the shown case the airfoil initially encounters a steady flow at zero incidence 
(a- = 0 deg), then pitches by 5 deg, and finally holds this angle ( a f = 5 deg). The airfoil piv- 
Ots on the leading edge, and the pitch-up motion takes the same time the airfoil needs to 
travel one chordlength (tv- 1). A Reynolds number of 1 million was chosen and the tran- 
sition constant G y was set equal to 120. To ensure a proper resolution of the spatial and 
temporal scales tife following gridsizes were selected to be in effect for the pilot runs. 1 he 
external inviscid flow is determined at 100 points distributed along the airfoil contour 
Boundarv laver profiles are computed at 90 streamwise stations on each of the upper and 
the lower surface. Each profile consists of 35 to 75 velocity vectors, depending on the 
boundarv laver thickness. The grid in normal direction is self-adapting in that the code au- 
tomatically figures whether to increase or to decrease the number of points across the lay er. 
Both potential and boundary layer flows are determined at 200 discrete time-levels, of which 
the first 100 are distributed over the ramp, and the second 100 over the stationary part of 
the airfoil motion. The time step is further controlled by a CFL-hke-condition, which, in 
case of a flow reversal, might enforce shorter time steps than chosen by the program. This 
constitutes roughly the input, which is provided to the program for each single case. 1 here 
is other data, like tolerance and relaxation parameter, which is set in the code. 

The graphic output of the program includes snapshots of the boundary layer velocity pro- 
files and the pressure distributions at selected time-levels, carpet plots for the distributions 
of skin friction and integral boundary layer thicknesses, and a time history of characteristic 
points like those of zero skin friction and the onset of transition. Figure 1 and its contin- 
uations illustrate the change in the velocity profiles and in the pressure distribution over the 
period the airfoil travels 5 chordlengths. The shown plots refer to the time instances at 
which the airfoil has covered a distance of 0.0, 0.5, 1.0, 1.5, 2.0, 3.0, 4.0, and o.O 
chordlengths, respectively. This is indicated by a marker in the angle of attack history at 
the lower left comer of each plot. The upper parts of the plots show the streamwise v eloc- 
ities inside the boundary layer at various stations on the airfoil contour. Additionally, the 
locations of the begin and end of transition are marked by a circle to allow for an easy 
reading of the different flow regimes. Since the boundary layer thickness amounts to few 
percent of the chord at the most, the vertical scales of the velocity profiles are enlarged by 
a factor of 10. The pitch-up affects the boundary layer primarily in two ways: first, the 
transition on the upper surface is shifted toward the leading edge, whereas on the lower 
surface toward the trailing edge. Second, the boundary layer grows on the upper surface, 
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in particular in the aft portion of the airfoil. Corresponding to the shift of transition, the 
regions of laminar and turbulent flow get redistributed: The flow starts out with 40 percent 
laminar, 15 percent transitional, and 45 percent turbulent flow on both surfaces. When the 
flow settles again after the pitch-up, the upper surface shows 10 percent laminar, 5 percent 
transitional, and 85 percent turbulent flow, while the lower surface shows 80, 10, and 10 
percent. Finally, there is also a good correlation between the tvpe of flow and the shape 
of the velocity profiles (almost textbook-like). The turbulent profiles can be identified by 
the rapid flow acceleration taking place very close to the wall. In contrast, the laminar 
profiles exhibit a much slower, yet more persistent, flow acceleration. The lower parts of 
Figure 1 present the pressure distribution and the location of the stagnation point with re- 
spect to the airfoil contour. Suction is indicated by arrows pointing from the airfoil, pres- 
sure by arrows pointing to the airfoil. In the course of the pitch-up the suction on the'upper 
surface becomes much stronger, the initially symmetric pressure around the stagnation 
point shifts to the lower surface, and the suction on the lower surface reduces to the static 
pressure of the freestream. The peculiar pressure rise just upstream of the trailing edae is 
due to the neglection of the displacement effect in the inviscid flow calculation. Both the 
results for the boundary layer and the inviscid pressure field considerably lag behind the 
motion. When the airfoil first reaches its final angle of attack, the boundary layer profiles 
and the pressure distribution are no way near their final, stationary values. Even after four 
extra traveled chordlengths the results are not perfectly steady/ though they then have 
pretty much approached their final state. 


The following series of five plots consists of carpet plots for the external velocity distrib- 
utions, the distributions of skin friction, displacement thickness, and momentum thickness, 
and the intermittency distributions. The development of the inviscid velocities along the 
airfoil contour with respect to time is shown in Figure 2. What has been said for the pres- 
sure applies equally to the velocities. They are only reproduced here, because thev repre- 
sent the input for the boundary layer flow solver. The most important output of a boundary 
layer calculation is (for practical purposes) the distribution of skin friction. Onlv this in- 
formation allows the evaluation of the viscous drag of a body. The coefficient of/kin fric- 
tion c,- is a measure of the viscous shearing stress at the wall in terms of the dynamic 
pressure and can be computed from the velocity profiles according to 


c f = 
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The variation of this coefficient with the chordwise distance and as a function of time is 
depicted in Figure 3. The key message of this figure consists in that the locus of transition 
exercises a predominant influence on the viscous flow results. More than everything else, 
transition shapes the development of skin friction. The steep increase in skin friction, which 
shows up in virtually every single curve of Figure 3, is associated with transition. This be- 
havior is entailed by the higher velocities of fluid particles near the wall in turbulent 
boundary layers, as compared to laminar ones. The shifting of transition due to the 
pitch-up is reflected as well in the shifting of the ramp-type increase in skin friction. The 
cur\es for skin friction (and the other boundary layer results) do not go all the wav up to 
the trailing edge. Stations, whose boundary layer thickness had exceeded a certain limit, had 
to be dropped from the calculation for stability reasons. Convergence problems, as experi- 
enced for the last two stations in the steady calculations, are attributed to the direct mode 
of the boundary' layer code, which cannot cope with flow separations. It should be pointed 
out, however, that these flow separations only occur because of the neglection of the dis- 
placement effect in the inviscid calculation. The inclusion of this effect shall take care of the 
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unrealistic strong pressure rise just upstream of the trailing edge and consequently avoid 
this kind of flow separation. The rapid growth of the boundary layer at the most down- 
stream station is believed to be related to the incomplete posing of the problem in that no 
boundary conditions are prescribed over the reversed flow region at the downstream face 
of the computational domain. 

That leads us to the ways of measuring the thickness of a boundary layer. Not only be- 
cause of the ambiguity involved in the definition of a "linear boundary layer thickne s 
integral thicknesses based, for example on the deficit of mass flow ™°^ UI ^ h du u e s ^ 
the "presence of the boundary layer, have been introduced. The most commonly used 
measure for the boundary layer thickness is the so-called displacement thickness d , whic 

is defined as 


il 
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(4.3) 


The displacement thickness can be interpreted as the p-S 4 

streamlines are shifted outwards owing to the presence of the boundary layer, rigure -4 
shows the variations of this thickness as a function of the chordwise distance and time. 
Transition also leaves visible marks on the development of the displacement thickness in 
that it is accompanied bv a short cessation of the otherwise continuous growth of the dis- 
placement thickness. This is caused by the sudden entrainment of faster fluid particles mto 
the regions close bv the wall, associated with the first occurrence of turbulence. Not sur 
trend canno’t hold very long, and particularly, 
bulent boundary layer begin to spread further into the initially inviscid flowfleld. This can 
he seen in the resumed now even stronger growth of the turbulent boundary layer shor 
downstream of transition. While the uneven displacement of the flows on the upper and the 
lower surface is responsible for the loss of lift, as suffered in viscous flow compared to the 
S another integral thickness that relates to the drag, exerted on a body 

m viscous flows In a similar way the displacement thickness is based on the loss of mass 
flow, the loss of momentum in a boundary layer can be used to define a second Integra 
thickness, the so-called momentum thickness 6 
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The development of the momentum thickness distributions in time is illustrated in Figure 
5 The historically substantiated importance of the momentum thickness has been reduced 
Ince the ntroduct on of finite difference methods, because flow reversals need no longer 
* mdfcttd bv the shape factor. The last of the carpet plots refers to the development of 
the taermuency distributions, which is shown in Figure 6. The intemuttency distribution 
acco'unTsTor^hJ transitional flow region, which exists between the laminar and the fully 
turbulent flow regimes! The most striking efTect in Figure 6 is the steepening ; of the in ter- 
mittencv distribution, when transition moves toward the leading edge, and the flattem g, 
whenTransition moves toward the trailing edge. This corresponds to a shortening or a 
lengthening of the transitional flow region, respectively. Such a behavior is not surprising 
in view of The role the external velocity plays in the equal, on lor the intermmency. 

picture 7 summarizes the flow results in terms of time histories for characteristic points. 
The presented information includes the points of zero skin friction 

and the stagnation point. The reader is reminded of two earlier stated findings. First, with 
the exception of steadv flow the point of zero skin friction does not coincide with the pom 
tf tow ^s P eparaln and second, the occurrence of flow reversals is partly attributed to the 
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neglection of the displacement effect. Even though, the historv of the zero skin friction 
points shows a rather interesting behavior Recalling that the airfoil pivots on the leading 
edge, the downward motion of the aft portion of the airfoil should give rise to a decreased 
entrainment on the upper surface and an increased entrainment on the lower surface. This, 
m turn, results in an extension of the reversed flow region on the upper surface and in the 
disappearance of it on the lower surface during the pitch-up. Most remarkable is the be- 
havior when the airfoil stops its rotation. The flow, that has gotten used to the effects of 
pitching, reacts to the loss of these effects correspondingly, namely with the reduction of 
the reversed flow region on the upper surface and the development of a flow reversal on the 
lower surface. Thereafter the flow slowly approaches the steady results with no flow reversal 
on the lower and 4 percent on the upper surface. As far as transition is concerned, the re- 
sults for upper and lower surface primarily differ in the amounts they lag the airfoil motion. 

1 he transition point on the upper surface reaches its steady-state position way sooner than 
the transition point on the lower surface does. This is related to the resistance an originally 
turbulent boundary layer offers to the transitioning to the laminar regime. 

Though the results demonstrate that the method can qualitatively predict the essential flow 
features, they do emphasize the need for including an interaction model. Further, in order 
to quantitatively validate the code, a comprehensive comparison with an experimental in- 
vestigation, that provides extensive data of the boundary layer, will be necessarv. 
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5. Concluding remarks 


This report describes the first step towards the development of a 1^''”““'! 

method for unsteady flows. Such methods require the successive solution of the ext ^nal 
inviscid flow and of the viscous flow in the boundary layers, linked together by a suitable 
interaction model. Our work has progressed towards the completion of the inviscid and 
viscous components of the method such that in the present status the code can be run in 
the so-called direct mode, which involves the solution of the viscous and the invisc d flo s 
according to the classical, uninteracted conception. The chosen methodology utilizes a 
unsteady panel method for the solution of the inviscid part of the flow, and a finite defer- 
ence method for the viscous part of the flow. The panel method derives from the stead 
technique of Smith and Hess, but additionally accounts for the continuous sheddin 0 of 
vorticitv into the trailing wake. The boundary layer method adopts Cebeci s technique of 
regular and characteristic boxes, which allows a stable integration of the unsteady^bounda^ 
laver equations in presence of flow reversals. A large portion of this report is devoted to the 
de'scriptionof these techniques, whereas only little room is given to the presentation and 
analysis of results. The findings of the limited testing of the code and the lessons learned 
during the development of the code can be summarized in the following, yet preliminary 

conclusions: 

1 In view of the influence the locus of transition exercises on the viscous flow results the 
modeling of transition becomes the key factor in the prediction of unsteady nscous 
flows. It seems at least doubtful, that a simplified model, like the currently used Mic 
criterion, can cope with such a challenge. On top of it, the numerical implementation 
of the steadv Michel criterion proved difficult in the unsteady environment. 

2 In case of an open region of reversed flow, the boundary layer exhibited rapid growth 

at the most downstream station. It is believed that this behavior is related to the in- 
complete posing of the problem. The presence of reversed flow over part of the 
downstream face of the computational domain necessitates the prescription of bound- 
arv conditions, corresponding to the upstream boundary' conditions at the stagnation 
noint there as well. Such conditions are not known a priori. , 

3 The occurrence of streamwise oscillations in the near-wall solutions of the boundary 
laver demanded a modification of the original box method. Since the averaging, re- 
sulting from the centered discretization, was identified as the cause of the problem the 
scheme was modified in that the solution is now taken midway between succeeding, 

4 The a SSqufoTregular and characteristic boxes is believed to be faulty with regard 
to the dSe e rufe. The damage caused by not accounting for the proper domain 
of dependence can hardly be evaluated, if only because of the lack of a faukless refer- 
ence scheme. The flawless implementation of the dependence principle would likely re- 
quire Se abandonment of the characteristic box scheme n view of the unargued 
merits of this scheme, such a step requires particularly careful consideration. 

5. Sufficient evidence has been found to endorse the argument that in tajers *c 
concept of stability is separate from the dependence principle. Th s is not in contra 
diction with the original CFL-condition /15/, which does imply such an identity since 
the CFL-condition strictly applies to linear, hyperbolic equations only, which are 

solved bv means of an explicit scheme. . . . . 

6 Considering the magnitudes of the boundary layer thicknesses witnessed in the calcu- 
lations thus far, the assumption, that the normal pressure gradient is negligible m the 
boundary layer, might not hold. The removal of such an assumption, however, would 
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mean to sacrifice the computational efficiency of the space-marching boundarv laver 
approach. 

A conceptual comparison with related studies does not reveal anything fundamentally new. 
The investigations of most relevance to our approach are those by Henkes and Veldman 
(on the breakdown of unsteady interacting boundary layer description / 23/), and by Geissler 
(on unsteady separation under dynamic stall conditions / 20/). With regard to the above 
addressed criticism, those two investigations do not do a better job either. Both “cheat" on 
transition: The former study deals with laminar flows only, the latter one assumes the onset 
of transition right at the stagnation point. Neither of the studies provides answers to the 
questions of open separation and the fulfillment of the dependence rule. In order to com- 
pute flow reversals Geissler adopts the zig-zag scheme, which, from a physical point of 
view, in our opinion is less favorable than the characteristic box scheme.' For the same 
purpose Henkes and Veldman employ the Crank-Nicolson method with upwind- 
differencing of the streamwise convection term. Concerning the interactive coupling of the 
viscous and inviscid components of the flow, however, both Henkes and Veldman, and 
Geissler are one step ahead of us, since both have such a feature already included in their 
present codes. That brings us to the improvements, which we anticipate to achieve in the 
short term, and which, in a broader context, we believe are necessarv to advance the current 
abilities of predicting viscous flow’s: 

1. Our method will now be extended to include a full interaction model, which should then 
permit the calculation of (moderately) separated flows. This goal shall be accomplished 
in two steps. First, the potential flow solver must be modified to account for the dis- 
placement effect due to the presence of the boundary layer. Second, the boundary laver 
method must be adapted to allow for the unknown status of the external velocity' in the 
calculation of the boundary layer profiles. 

2. The second near-term goal of our efforts refers to an improved implementation of the 
dependence rule, as well as of the law of forbidden signals. Though such an endeavour 
might not be beneficial from a practical point of view, w r e consider it important to rec- 
oncile the numerics with the underlying physics. 

3. The extension of the interactive boundary 'layer description to the wake will not only 
resolve the problems with the unknowm, yet (in case of a flow reversal) necessarv', 
downstream boundary conditions, but will also greatly advance the ability to handle 
separated flows. 

4. There is a great need to improve the modeling of transition and turbulence. In the long 
term, only a thorough understanding of these two phenomena will lead to a more reli- 
able prediction of viscous flows. 
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Appendix A. The regular box scheme 


This appendix provides detailed information on the discretization of the unsteady boundary 
layer equations by means of the regular box scheme. With the exception of reversed flow 
regions and the immediate neighborhood of the stagnation point, the regular box scheme 
is applied throughout the fiowfield. The box method is a midpoint scheme, which expresses 
all functions or derivatives in terms of quantities at the corners of a computational cell. To 
this end the boundary layer equations are written as a system of three first-order partial 
differential equations ‘ 


-kw + FV- 


8F_ _ dU_ 
drj dt, 

dU TJ 8JU 
dr 3S, 


= 0 


HL_ W A1L 

dr d{ 


(A. 1) 


(A.2) 


This system is complemented by the usual boundary conditions, which are the no slip 
condition at the wall and the prescription of the external velocity at the boundarv laver 
edge 
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V = Ve- V=W{t,r) (^- 4 ) 


For the finite difference approximation a rectangular mesh with computational cells of size 
r, x hjX k n is introduced. The position in space and time is denoted by a treble suffix, for 
example t£", where; is referring to the cross-stream direction, i to the streamwise direction, 
and n to the time. The equations can now be approximated in terms of central differences 
and two-point averages, involving the values at the comers of the computational cell only. 
Accordingly, the continuity equation (AA) is evaluated at the point 2 . rj t ) the 
momentum equation {A.2) at the center of the cell 2 , n HU2 , and t£e auxiliary 

equation (A. 3) at the point (?,.^_ l/2 , rj 
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The boundary’ conditions read in terms of nodal values according to 

U l { n = 0, F\' n = 0 

U l j n = W‘' n (A. Si) 
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Since the resulting equations are algebraically nonlinear in the unknowns LP/ , V)’ , F/ , the 
scheme proceeds in linearizing the equations by Newton's method. This procedure ap- 
proaches the solution in an iterative process, in which the variables of the kiH step are set 
equal the value of the previous iteration (k-1) plus a small correction 

Uj' n,K = U l j’ n ' K ~ ] + SUj-’ n ' K where S U‘' n ’ K 4 U)' n ' K 1 

yi,n,K _ yi/i,K—i _|_ where S vj’ n ’ K Vj’ n ’ K 1 (A. 9) 

pl.n, k = pi.n. K-1 + s fj.n. k where rfA.* ^ «-> 

After dropping the terms involving quadratic expressions of the correction b and collecting 
the like-terms, the unsteady boundary layer equations become 
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This final set of algebraic equations is Unear in the unknown corrections dbj , 
SV i ' n ' K , SFj' n ’ K . To avoid lengthy expressions equations (T.10) and (/4.11) make use of the 

following abbreviations 
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Together with the boundary conditions equations (4.10) through (AA2) constitute a system 
of j/ algebraically linear equations (where J is the number of gridpoints across the bound- 
ary layer). The equations are arranged such that the resulting matrix is of the block- 
tridiagonal type. Each line of that matrix refers to a set of three equations (3 x 3 blocks) 

1. The first set consists of the following equations 

a. The boundary condition requiring zero tangential velocity at the wall 

b. The boundary condition requiring zero normal velocity at the wall 

c. The auxiliary equation (A. 12) 

2. Other than the first and last sets include the following 

a. The continuity equation (,4.10) 

b. The momentum equation (,4.11) 

c. The auxiliary equation (,4.12) 

3. The last set is made up by the following three equations 

a. The continuity equation (,4.10) 

b. The momentum equation (,4.11) 

c. The boundary condition prescribing the external velocity at the boundarv laver 
edge 


The solution of this system is obtained by the block elimination method, which consists of 
two sweeps. The forward sweep eliminates the blocks below' the main diagonal. In the 
backward sweep the solution is calculated from recursion formulas. The actually employed 
method is specifically designed for the above described system, such that it is more efficient 
than a general algorithm for block tridiagonal systems. 
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Appendix B. The characteristic box scheme 

This appendix deals with the difference approximation to the unsteady boundary layer 
equations bv means of the characteristic box scheme. The characteristic box scheme was 
introduced bv Cebeci and Stewartson in order to allow for some upstream influence when 
the local velocity is directed upstream. Consequently the characteristic box scheme is ap- 
Dlied in regions of backflow and in the immediate neighborhood of the stagnation point. 
Further the calculation of the upstream boundary conditions involves the iterative use of 
a modified characteristic box method, in which the continuity equation is decoupled from 
the momentum equation (see Appendix C.). In order to utilize the box scheme the un y 
boundary layer equations, with the streamwise convection term being evaluated along a 
subcharacteristic, are reduced to a system of first-order partial differential equations 
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Here the newlv introduced coordinate s denotes the distance along the subcharacteristic. 
Onlv the momentum equation (5.2) differs from that of the regular box scheme while the 
continuity equation (5.1), the auxiliary equation (5.3), and the boundary conditions are 
identical to those of the reeular box scheme. The discretization closely follows the proce- 
dureas given in Appendix A. The centering of the continuity and the auxiliary equation 
remains unchanged, namely the continuity equation is evaluated at the point 
(c -> n , T„), and the auxiliary equation at the point (<;,, r\ } . 1/2 , t„). The difference mole 
cule foMite momentum equation is here a rectangular plane, generated bj 'the, - -axis and 
the direction of the subcharacteristic. Hence the momentum equation is centered at the 
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The suffix m denotes the streamwise location where the subcharacteristic intersects the 
previous timeline. The station % p is located halfway between $,■ and Any ow qu y 
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which needs to be evaluated at <f m or <N, must be interpolated from the surrounding nodal 
values. The currently employed interpolation algorithm utilizes a Lagrangian polynomial. 

Together with the boundary conditions equations (5.4) through (5.6) constitute a svstem 
in the unknowns Ly, V'y, and Py. Because of the nonlinear convection term Newton's 
method is adopted to linearize the equations. The linearized equations can be written in 
terms of the unknown corrections 5 by*, and 5py K 
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where the superscript k denotes the iteration counter. The structure of the resulting svstem 
is identical to that of the regular box scheme, and the only change occurs in the coefficients 
which now are defined bv 
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Since the resulting equations of the characteristic box scheme are of the same form as those 
for the regular box, the same algorithm can be employed to obtain a solution. The 
equations are arranged in the same manner as described in Appendix A in order to form a 
block-tndiagonal system, which is solved by Keller's block elimination method. 
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Appendix C. The upstream boundary conditions 


As is well known a mathematically well posed problem requires that the governing differ- 
ential equations be accompanied by appropriate initial and boundary conditions. This ap- 
pendix describes the method which provides the upstream boundary conditions for the 
unsteady boundary layer equations. The proposed method consists of an iterative compu- 
tation of the flow at the stagnation point, and the adjacent upper and lower surface points. 
For the computation of the flow at the stagnation point, the method adopts a modified 
characteristic box scheme. In the following we will first derive the difference equations for 
the modified characteristic box scheme, and we will then outline the proceeding how to 
obtain the upstream boundary conditions. 

The modified characteristic box scheme differs from the original one in that the continuity 
and the momentum equations of the modified scheme are solved in a decoupled manner. 
Consequently, the tangential component u and the normal component o of the velocity are 
calculated separatelv, the first one by solving the momentum equation, the latter one by 
solving the continuity equation. It is this decoupling and the use of a characteristic box 
scheme, which allows" the solution of the flow at the very first point on a new timeline. 1 he 
difficulty of this endeavour lies in that there is absolutely no information available on the 
new timeline. This disables the regular box scheme, because it would involve a 
C - derivative on the new timeline (requiring a known point there). The solution ot the 
continuity equation is impeded for the same reason. However, there is a drawback associ- 
ated with the decoupled method, that is the normal velocity component o must be assumed 
known in the momentum equation. This entails an iterative procedure, which proceeds as 
follows. The first step of each of those iterations solves the momentum equation for a 
known though provisional value of the normal velocity component o, which is taken from 
the previous level of iteration. Subsequently, the second step updates the normal velocity 
component o bv solving the continuity equation. This procedure is repeated till successive 
values of the normal velocity component o differ by less than a prescribed small quantity. 

We will now briefly present the equations for the modified characteristic box scheme. The 
system of first-order partial differential equations is reduced to two (since the continuity 
equation is left out), namely the momentum equation and the auxiliary equation, which 
results from the introduction of the new dependent variable V 
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ML _ v = o (C.2) 

8n 

The tilde indicates provisional values that are taken from the previous level of iteration. 
Because F is treated as a known coefficient, there are only two boundary conditions neces- 
sary 


n = 0: U - 0 

n = n; u~W{t,x) 


(C.3) 


The discretization follows the previously outlined procedure, with the momentum equation 
being centered at the point (c p , anci ^ auxiliary equation at the point 

(£ ;) rj _ u2 , t„). After employing Newton's linearization the system can be rewritten as 
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Together with the boundary conditions equations (C.4) and (C.5) form a system of 2J al- 
gebraically linear equations (where J is the number of gridpoints across the boundarv 
layer). After proper arrangement of the equations, the resulting block tridiagonal svstem 
can be solved by the block elimination method, which has here to operate on 2 x 2 blocks. 

We will finally turn to the procedure of obtaining the upstream boundarv conditions. The 
flow at the stagnation point, which represents the upstream boundary conditions for both 
upper and lower surface, is obtained by an iterative procedure, in which each sweep in- 
\ol\es the solution for the flow at the stagnation point and at the adjacent upper and lower 
surface points. In detail the procedure consists of the following steps: 

1. For the initial iteration F must be guessed at the stagnation point. For that purpose F 
is linearly extrapolated on basis of the two preceding timelines 

- ■'7' + i~- 7' - 7) (c.7) 

— I 

The superscript outside the brackets refers to the number of the sweep, the 0 in the 
above equation denotes the initial sweep. 

2. The flow at the stagnation point is solved by means of the modified characteristic box 
scheme. The results b and Fare of approximate nature, since Fis assumed known from 
the previous level of iteration. 

3. The flow at the adjacent upper surface point is solved by the usual procedure. That is 
the regular box sheme is used for that part of the boundary layer, where the tangential 
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velocity component is directed downstream, while the characteristic box scheme is em- 
ployed in regions of backflow. 

4. The flow at the adjacent lower surface point is solved by the same procedure that is 
used for the upper surface point. 

5. The provisional values of F are updated by utilizing the (so-far unused) continuity 
equation 


[?/ n y = [?at + 




(C.8) 


subject to the boundary condition (of zero normal velocity at the wall) 

[?{"]* = 0 (C.9) 

The £ — derivatives in equation (C.8) are approximated by means of Lagrangian 
polynomials, which are based upon the values at the stagnation point (index /), the 
upper surface point (index / + 1), and the lower surface point (index / - 1) 
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6. Finallv it is checked, whether the values of F have converged to a prescribed accuracy. 
If successive values of F at the boundary layer edge differ by less than a given tolerance, 
then the procedure is exited and a set of upstream boundary conditions has been es- 
tablished. Else the procedure continues with the second step and is repeated till the 
criterion of convergence is met. 
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The solution of the nonlinear boundary layer equations requires that the equations be 
linearized by Newton s method. In this method the dependent variables are substituted bv 
the (known) variables of the previous level of iteration plus a small correction. After drop- 
ping the quadratic terms, the equations are solved for the unknown corrections. In order 
to start the iteration an initial guess of the dependent variables must be provided. The 
quality of this initial guess has a great deal of influence on the number of the required it- 
erations, and ultimately, a bad guess can lead to diverging results. A simple, though inac- 
curate guess would be given by the values of either the adjacent upstream station or the 
previous timeline. In our code the initial guess for U is obtained from a two-step procedure, 
which consists of an extrapolation and an adjustment to the imposed external velocity. The 
extrapolation is preferably performed in time, since the results of two succeeding timelines 
are usually closer than those of two adjacent streamwise stations. The velocity profile ob- 
tained from the extrapolation is adjusted to the given external velocity in order to auto- 
matically satisfy the boundary condition at the boundary layer edge. Guesses for the 
remaining two variables V and F follow from the solution of the auxiliary and the conti- 
nuity equations with the guessed U profile as input. 

If possible we base our scheme on a linear extrapolation of U in time and a proportional 
adjustment of the extrapolated profile (designated by a caret) to the external velocity 
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This extrapolation scheme can only be employed, if both W’> n and IS/’ 0 do not equal zero. 
In case one of these conditions is violated, we adopt a slightly modified scheme. There the 
adjustment to the external velocity is achieved by the choice of a proper dUjdz 
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This scheme requires that both terms ( W’ n - W 1 *-') and - iv'' n ~ 2 ) be nonzero. If 

both procedures fail, we combine the first step of the first procedure with a single adjust- 
ment of the outermost point to the given external velocity. 

The values of V can be guessed by differentiating the above obtained U 
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An estimation of F is obtained from the substitution of the guessed U into the continuity 
equation (subject to the condition of zero normal velocity at the wall) 
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Appendix E. On instabilities of the near-wall solutions 

The standard box method, as devised by Keller and Cebeci, allows the occurrence of 
streamwise oscillations in the near-wall solutions of the boundary layer. The oscillations 
are caused by the centered discretization of the momentum equation, which, by taking 
two-point averages, may offset the errors of symmetrically overshooting solutions at the 
upstream and the downstream faces of the difference molecule. Similar observations were 
made by Drela and Thompkins / 17, 1 8/ in their study of non-unique solutions of boundary 
laver equations. Though the possible occurrence of these oscillations is, according to our 
experience, an inherent shortcoming of the box method, Cebeci never reported on problems 
of this kind. This might be related to the fact that a well-smoothed external velocity dis- 
tribution (as input) and the use of a well-adapted boundary' layer grid do suppress the os- 
cillations to a areat extent. As a matter of fact, we only became aware of the problem, when 
we tried to run a steady flow case in physical coordinates. In the following we will give the 
reason for the occurrence of overshooting solutions, and biefly describe our proposed 
modification to the box scheme in order to ensure a smooth solution. 

For the sake of simplicity we shall explain the development of these instabilities for the 
steady, laminar boundary layer. The momentum equation reduces then at the wall to 


-.2 

C U 

dy 2 


S - u a 


ou„ 


ox 


(£. 1 ) 


This equation retains only the viscous and the pressure terms, and is exact at y — 0. Further, 
the equation approximately holds true in the near-wall field, where those two terms domi- 
nate the balance of momentum. The box methods requires the momentum equation be 
discretized at the center of the box such that equation (£.1) is approximated by 
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Since this equation only deals with the average of d 2 u/cy ■ between two successive 
streamwise stations, d 2 uldy 2 , and consequently u, may exhibit oscillatory- deviations and still 
satisfv the equation. The error at the ith station is simply compensated by an equal error 
of opposite sign at the i-lst station. Unfortunately, at the wall there is no mechanism to 
dampen these oscillations, whereas there is one further away from the wall. With the 
convective terms gaining weight in the momentum balance, the profile of u is increasingly 
constrained by the need to approach the prescribed, external velocity. Finally there remains 
the question, how the oscillations get started. Two causes are believed to initiate the first 
deviation from a smooth solution: First, a slightly uneven external velocity distribution, as 
was demonstrated by Drela with his zero-pressure gradient flow responding to a o percent 
jump of the edge velocity, and second, the sudden change of How conditions at successive 
streamwise stations on t] = constant lines, as possible in poorly adaptiv grids that do not 
account for the growth of the boundary layer. Before we proceed with our proposed mod- 
ification of the box scheme, we would like to point out again that these instabilities do only 
occur near the wall, and they have not been observed for the standard, steady scheme with 
the use of transformed coordinates (at least the amplitudes of the oscillations there are re- 
duced to a negligible degree). 

The oscillations are numerically only possible, because the point, where the boundary- layer 
equations are satisfied, is separate from the point, where the solution is actually computed. 
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Based on this finding the box method was modified in that the solution is now taken at that 
point, which actually satisfies the equations, which is the center of the box. The box method 
still solves for the variables at the downstream face of the difference molecule, but the 
variables, which are given as the solution, and which are used in the further calculations, 
are those associated with the center of the difference molecule. These variables are the 
two-point averages of the upstream and the downstream faces of the difference molecule in 
case of steady flow, or the four-point averages of the upstream and the downstream, as well 
as the current and the previous timeline nodes in case of unsteady flow. Though this pro- 
cedure completely eliminates the occurrence of streamwise oscillations in the near-wall sol- 
utions, it ma> be a mixed blessing. While the modified scheme formally retains second order 
accuracy, the actual truncation error might become larger due to the implicit smoothing. 
Furthermore, under the assumption of the same meshsize in the streamwise direction, the 
required computational effort will double. In view of the shortcomings of both methods, 
one might well ventilate the question, whether a non-centered scheme will do a better job' 
Non-centered schemes will not experience instabilities of the above described kind, but in 
order to maintain second order accuracy, the streamwise derivatives must be approximated 
by one-sided, three-point difference representations. Veldman /45 / and his collaborators 
/23/ have made extensive use of the non-centered Crank-Nicolson method, and in view of 
their success a consideration of such a scheme might well be worthwhile. 
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BOUNDARY LAYER VELOCITY PROFILES 

Symbols: T1 ... Bogin of Transition 
T2 ... End of Transition 
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POTENTIAL FLOW PRESSURE DISTRIBUTION 
Symbols: SP ... Stagnation Point 
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Pivot: (0.00,0.0) 
Rise Time: 1.00 


Figure 1. The boundary layer velocity profiles and the pressure distributions of a NACA 
0012, undergoing a ramp-change in angle of attack. 


Results of a sample calculation 
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BOUNDARY LAYER VELOCITY PROFILES 

Symbols: T1 ... Begin of Transition 
T2 ... End of Transition 


B.L SCALES 





PRESSURE SCALE 


TIME HISTORY OF THE ANGLE OF ATTACK 
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NACA 0012 

87=11=10 

Reynolds No: 1 mill, 

Trans Const: 120., 

Bl Stations: 181 
Timesteps: 200 

RAMP CHANGE IN AOA 

Initial AoA: Odeg, Pivot: (0.00,0.0) 

Final AoA: 5deg, Rise Time: 1.00 


Figure 1. (continued) 
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Results of a sample calculation 




BOUNDARY LAYER VELOCITY PROFILES 

Symbols: T1 ... Begin of Transition 
T2 ... End of Transition 
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Trans Const: 120 , 
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Figure 1. (continued) 


Results of a sample calculation 
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BOUNDARY LAYER VELOCITY PROFILES 

Symbols: T1 ... Begin of Transition 
T2 ... End of Transition 
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Figure 1. (continued) 


Results of a sample calculation 








BOUNDARY LAYER VELOCITY PROFILES 

Symbols: T! ... Begin of Transition 
T2 End of Transition 




a 



NACA 0012 

87-11-10 


81 Stations: 181 
Tinri 0 step 8 : 200 

RAMP CHANGE IN AOA 

Initial AoA: Odeg, Pivot: (0.00,0.0) 

final AoA: 5daa. Rise Time: 1 00 


Figure 1. (continued) 
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Results of a sample calculation 
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Figure 5. The distributions of momentum thickness of a NACA 0012, undergoing a 
ramp-change in angle of attack. b 
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Results of a sample calculation 



Figure 6. 



Intermittency Factor, 7 


The distributions of intermittency of a NACA 0012, undergoing a ramp- 
change in angle of attack. 
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Results of a sample calculation 
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Figure 7. The histories of characteristic points in the flow past a NACA 0012, under- 
going a ramp-change in angle of attack. 


Results of a sample calculation 









